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Although quartz (a-form) is a mineral used in numerous applications wherein radiation exposure 
is an issue, the nature of the atomistic defects formed during radiation-induced damage have not 
been fully clarified. Especially, the extent of oxygen vacancy formation is still debated, which is 
an issue of primary importance as optical techniques based on charged oxygen vacancies have been 
utilized to assess the level of radiation damage in quartz. In this paper, molecular dynamics (MD) 
simulations are applied to study the effects of ballistic impacts on the atomic network of quartz. 
We show that the defects that are formed mainly consist of over-coordinated Si and O, as well as 
Si-O connectivity defects, e.g., small Si-O rings and edge-sharing Si tetrahedra. Oxygen vacancies, 
on the contrary, are found in relatively low abundance, suggesting that characterizations based on 
E' centers do not adequately capture radiation-induced structural damage in quartz. Finally, we 
evaluate the dependence on the incident energy, of the amount of each type of the point defects 
formed, and quantify unambiguously the threshold displacement energies for both O and Si atoms. 
These results provide a comprehensive basis to assess the nature and extent of radiation damage in 
quartz. 


I. INTRODUCTION 

Quartz (the dominant crystalline form of Si02 under 
ambient conditions) is one of the most abundant minerals 
in Earth’s crust. It is widely used in engineering appli¬ 
cations, e.g., in construction materials, semiconductors, 
integrated optics, ultrasonic devices, etc. Many 

of these applications require the material to withstand 
exposure to radiations. For example, when used as a 
concrete aggregate, quartz can be subjected to neutron 
radiation in nuclear power plants 017 ]. Therefore, it is 
important to establish a good understanding of the na¬ 
ture of the damages induced by radiations, so that the 
material durability and reliability can be improved. In 
addition, improved knowledge of radiation-induced dam¬ 
age mechanisms in quartz could also inform radiation 
damage sensitivities in similar materials such as silicate 
glasses, which are also of substantial interest for radiation 
related applications including nuclear waste encasement 
and optics in fusion reactors 0. 

Numerous experimental studies have been conducted 
to assess the nature and extent of damage induced in 
quartz by various types of irradiations, and often using 
neutrons and different types of ions [TJ 0 IMT^. It is 
generally admitted that different types of irradiation re¬ 
sult in similar damage to the atomic structure [12] . Such 
damage has been described in term of an accumulation of 
point defects, which eventually causes structural amor- 
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phization mm- However, the final amorphized struc¬ 
ture is found to be different from that of glassy silica, and 
the detailed nature of the formed point defects has not 
been clearly identified nnjE]. Numerous studies have 
focused on the oxygen vacancies created during irradia¬ 
tion IDITIHIS]. Indeed, it is well-known that the oxygen 
vacancies present in quartz after ionization or electronic 
radiations can trap electrons, forming the so-called E' 
center, a color center that modifies the optical proper¬ 
ties of quartz HZ]. This has been used to assess, by 
optical techniques, the level of damage and the thresh¬ 
old displacement energy (Ed) [Dill], a critical quantity 
commonly used to estimate radiation damage [18] . How¬ 
ever, the reported Ed is substantially larger for oxygen 
than for silicon. This appears counterintuitive just by 
considering that four bonds need to be broken to dis¬ 
place a silicon atom, but the rupture of only two of such 
bonds are needed for oxygen displacement [2]. The re¬ 
ported Ed for oxygen in quartz is also almost twice as 
large compared to that in fused silica HU, which is unex¬ 
pected. Additional experimental evidence suggests that 
E' center creation is accompanied by substantial forma¬ 
tion of other unidentified defects (up to 500 other defects 
for every E' center), suggesting that relying only on the 
oxygen vacancies to assess the damage and determine rel¬ 
evant properties may not be appropriate [T0|T2]. These 
complexities in structural damage models impede better 
understanding of the damage process, and also create dif¬ 
ficulties in accurately predicting the damage evolved in 
materials, for engineering purposes. 

Molecular dynamics (MD) simulations provide a 
unique ability to access the damage process in irradi- 
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ated quartz. In the past decade, several studies have 
applied MD to explicitly simulate radiation damage in 
silicates |19lE5] . For crystalline quartz, however, such 
studies have been limited. To this end, herein, we per¬ 
form MD simulations of radiation-induced atomic struc¬ 
tural damage in quartz. By detailed analysis of MD tra¬ 
jectories and careful consideration of the creation and 
annihilation oxygen vacancies, we confirm that oxygen 
vacancies only account for a small percentage of all the 
induced defects. The nature and quantity of the domi¬ 
nant defects and their dependence on the incident energy 
are evaluated. This information, together with improved 
estimations of threshold displacement energies for both 
oxygen and silicon atoms, allow for more accurate assess¬ 
ment of the radiation damage in quartz. Finally, we also 
examine damage to the connectivity of the Si-0 network, 
such as the modification of the ring size and the appear¬ 
ance of edge-sharing Si tetrahedra, which could serve as a 
structural basis for the radiation-induced amorphization 
at higher dosage. 


II. SIMULATION METHODS 

Molecular dynamics simulations were carried out with 
the LAMMPS code |26]. To ensure that both pristine 
and amorphized quartz can be properly simulated, we se¬ 
lected the BKS potential, m a widely adopted potential 
set for silica. Note that this potential has been previously 
used to study radiation damage in glassy silica and has 
been shown to provide realistic predictions HJ [531 1^ - 
The potential, although reliable at equilibrium separa¬ 
tion, does not provide any physical description of interac¬ 
tions when the atoms reach closer distances with respect 
to each other, as observed during high-energy collisions. 
These interactions can, instead, be described with the 
ZBL potential, which has been designed to tackle such 
atomic collisions [T3|. For a smooth transition, the BKS 
and ZBL potentials are linked by a polynomial function. 
The first and second derivatives are examined carefully 
when adjusting the polynomial transition to ensure that 
no additional extrema (that would result in unrealistic 
metastable states) are created by the process. The cutoff 
for the BKS potentials is set to 5.5 A, following the previ¬ 
ously reported protocol, [281129] and the ZBL potentials 
take effect below 0.8 A for Si-0, 1 A for 0-0, and 1 A 
for Si-Si interactions, respectively. 

Radiation damage is simulated within a supercell of 
quartz. The deposition of irradiation energy is simulated 
such that an atom in the supercell is randomly selected 
and accelerated with a given kinetic energy, correspond¬ 
ing to the impact (incident) energy. It must be noticed 
that the absorption of irradiation energy is not explic¬ 
itly simulated here. However, the kinetic energy used to 
accelerate the atom can be readily translated from the 
incident irradiation energy for ion implantation, or the 
incident radiation energy for neutron, electron or elec¬ 
trostatic radiations, if the absorption behavior is known. 


Therefore, for convenience, this acceleration energy will 
be referred to as the "incident energy" in the following 
sections. 

Following the energy deposition, the kinetic atom leads 
to a cascade of collisions with the neighboring atoms. A 
spherical region is created around the region that would 
be potentially damaged. Within this region, the tem¬ 
perature of the atoms is not controlled. A thermostat is 
instead placed on the rest of the system, outside of the 
spherical region, and set to 300K. Hence, the thermostat 
does not influence the true dynamics of the displaced 
atoms during the ballistic cascade, whereas the rest of 
the system, which is not directly affected by the cas¬ 
cade, serves as heat sink, ensuring that local heating and 
quenching effects are properly reproduced. Note that, 
since the damage region is dynamically created following 
the irradiation energy deposition, the need to maintain a 
specific area for heat sink is eliminated, allowing damage 
to be initiated from any atom in the system. 

During the initial cascade, a variable time step is in¬ 
voked to avoid excessive displacements of atoms during a 
single step, with a maximum time step of 1 fs. Depend¬ 
ing on the incident energy, most collisions occur within 
the first 1 to 3 ps, before the temperature of the system 
eventually stabilizes to 300K. Atomic trajectories during 
the first 3 ps are therefore analyzed in high-resolution to 
track the collisions. Afterwards, the simulation is con¬ 
tinued for another 50 ps, allowing the structure to relax. 
To obtain meaningful statistical data, this process is re¬ 
peated 100 times for each incident energy. At each cycle, 
an impact atom is randomly selected. Both the size of 
the supercell and the spherical region are selected based 
on the incident energy. For 300, 600 and 900 eV incident 
energies, 8x9x7 (4536 atoms), 10x12x9 (9720 atoms), 
13x15x12 (21060 atoms) quartz supercells, and 16 A, 
25 A, and 33 A radius spherical regions are used, respec¬ 
tively (see further details in the supporting information). 


III. RESULTS AND DISCUSSIONS 
A. Structural damage 

Radiation damage in oxides is commonly described by 
point defect formation miisni. For quartz too, this has 
been proposed as the primary damage manifestation at 
low dosage [TJ |TT] . Because the oxygen vacancy in quartz 
can trap localized charge, forming the E' (color) center, 
such defects have been used to study damage in neutron 
radiated and ion-implanted quartz [nini. However, we 
find that oxygen vacancy formation is not the main form 
of damage involved. 

In a cascade event (an example is shown in Fig. 
an atom absorbs the radiation energy and is accelerated 
(known as the Primary Knock-on Atom, or PKA). The 
high kinetic energy causes the PKA to substantially dis¬ 
place from its original position and, along the way, to 
collide with other atoms (Secondary Knock-on Atom, or 
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FIG. 1. An example 600 eV cascade trajectory. Silicon and oxygen atoms with a displacement greater than 1 A are shown as 
yellow and red spheres, respectively. Major high-velocity atoms during the cascade, such as the primary knock-on atom (PKA) 
and secondary knock-on atom (SKA), are highlighted with exaggerated sizes, and their initial positions are shown with beach 
balls. 


SKA). In the case of an oxygen PKA, such a collisions do 
indeed create a vacancy at its original lattice site, imme¬ 
diately after displacement, leaving two three-coordinated 
Si atoms (Si^^^). However, we observe that a significant 
number of these Si^^^ disappear within a fairly short time 
(<10 ps). The analysis of 100 cascades for the three se¬ 
lected incident energies reveals that, on average, less than 
7% of the PKAs eventually generate oxygen vacancies 
(and Si^^^ species) at their original positions, as shown in 
Fig. a). In addition, for all the three incident energies, 
even in the case of Si^^^ formation, no oxygen PKA is 
observed to simultaneously generate and maintain both 
of the two Si^^^ species that are predicted by the ideal 
senario. The same feature is observed when tracking all 
the knock-on atoms (KA) with a displacement greater 
than 4 A (note that the inter-tetrahedron distance in a- 
quartz is 3.06 A), as shown in Fig. ilb). 

Further analyses of the cascade trajectories reveal 
three mechanisms for the annihilation of the oxygen va¬ 
cancies following their initial creation. 

1. The first mechanism manifests by the recombina¬ 
tion of the interstitial and vacancy. This does not 
necessarily require the knock-on atom to return to 
its original position, but rather it can be achieved 
by a chain of atoms hopping to their neighboring 
sites, as illustrated in Fig. |^a). In our simulations, 
this has only been observed to occur when the dis¬ 
placement of the knock-on atom is smaller than 10 
A. However, its occurrence at 300~900 eV incident 


energies is already rare so it should not be a dom¬ 
inant mechanism for oxygen vacancy annihilation 
under higher incident energies. 

2. The second mechanism, which is most common in 
our simulations, manifests by a local relaxation 
around the oxygen vacancy. After the initial for¬ 
mation of the two Si^^\ they can relax away from 
the vacancy and form small rings with the nearby 
silicon tetrahedra. This involves the formation of 
three member Si-0 rings, as well as three-fold coor¬ 
dinated oxygen atoms, or tricluster oxygen, which 
are found to exist in some silicate glasses m- An 
example is shown in Fig. |^b), in which two three- 
member Si-0 rings are formed by the two Si^^^ relax 
into an interstitial site in the quartz structure. 

3. The third mechanism, i.e., a replacement collision, 
involves a primary atom knocking out but at the 
same time replacing a SKA. This mechanism is 
found to be more common for Si KAs, as will be 
discussed later. 

Other than oxygen vacancies, we note that other types 
of defects appear to form more readily during the cas¬ 
cade. As shown in Fig. Sa), a Si PKA almost always 
creates one or more under-coordinated oxygen atoms or, 
in other words, non-bridging oxygen atoms. For Si SKA, 
such probability decreases significantly due to the in¬ 
creased chance of replacement collisions and recombina¬ 
tion with nearby interstitial defects, but the probability 
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(a) 


(b) 


FIG. 2. Probabilities of creating zero or one three-coordinated Si (Si^^^) by (a) a O PKA, or (b) by all the knock-on atoms 
with displacement greater than 4 A generated by an O PKA. In any case, no simultaneous creation of two Si^^^ is found. 


Chain hopping 


Relaxation 



FIG. 3. Mechanisms for oxygen vacancies removals, (a) shows a series of neighboring hops to achieve the recombination, (b) 
shows the relaxation of Si^^^ around the oxygen vacancy to form small Si-0 rings and regain four-fold coordination. 


of forming non-bridging oxygens is still greater than 50% 
for all three incident energies (see Fig. I^b)). At the 
end, the incident energy does not significantly impact the 
overall probability of creating dangling oxygen atoms by 
the knock-on atoms. It is worth noting that Si knock-on 
atoms can also lead to formations of small Si-0 rings. 

Based on the results above, we note that different in¬ 
cident energies create similar types of point defects, i.e., 
under- and over-coordinated Si and O atoms. As illus¬ 
trated in Fig. [2 a knock-on atom usually creates under¬ 
coordinated defects at its original position. On the other 
hand, at the end of the cascade, displaced atoms usually 
remain in interstitial positions and form over-coordinated 
defects, such as tri-cluster oxygen and five-coordinated 
Si atoms. For a typical cone-shaped damage volume, at 
the end of the cascade, more over-coordinated defects are 
formed at the end of the cascade than under-coordinated 
defects around the initial knock-on atom. Indeed, the 
enumeration of all types of defects formed in a single 


cascade confirms such conclusion. As shown in Fig. 
oxygen vacancies (same quantity as Si^^\ as a knock-on 
atom dose not creates more than one oxygen vacancy) 
only constitute a small portion of the total defects for all 
incident energies. Also, we note that the ratio between 
different defect species is not a linear function of the inci¬ 
dent energy. Indeed, the formation of over-coordination 
defects appears to be far more pronounced than that of 
the under-coordination ones. Therefore, the relative per¬ 
centage of the oxygen vacancies, among all the point de¬ 
fects, decreases with the incident energy. As such, the 
E' center (charged oxygen vacancy) concentration alone 
may not be a proper indicator of radiation damage in 
quartz as illustrated in Fig. ib) below. 


The above conclusions agree with experimental obser¬ 
vations based on the energy loss during ion implanta¬ 
tion (incident energy >30 keV), according to which up to 
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(a) (b) 


FIG. 4. Probabilities of creating zero, one, two, three or four under-coordinated O (O^) by (a) a Si PKA, or (b) by all the 
knock-on atoms with displacement greater than 4 A generated by a Si PKA. 



FIG. 5. Average number of point defects created during one cascade is shown in (a). Data are fitted by power functions Nd= 
AK™, represented by dashed lines, (b) shows the ratio between number of each defect type (Nd) and number of Si^^^ (Ngjiii), 
extrapolated to higher incident energies using the fitting results. 


500 other defects could be created for every E' generated 
duns]. However, this estimation may be imprecise as 
the nature of the defects, and their formation energies, 
was not precisely known. To provide a more accnrate es¬ 
timation, we establish an analytical relationship between 
the number of defects (Nd) and the incident energy (E), 
which can then be nsed to extrapolate the data displayed 


in Fig. I^to higher incident energies. Althongh the lim¬ 
ited range of system sizes makes it difficult to determine 
the exact fnnctional form of the relationship, it appears 
that they can be fitted by a power fnnction of the form 
Nd= AE'^. The corresponding parameters are listed in 
TABLE H] Hence, the calcnlated fraction of oxygen va¬ 
cancies, with respect to the total amonnt of defects, de- 
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TABLE I. Fitting parameters for the energy dependence of 
the number of defects shown in Fig. |^a) (see text). 


A (IQ-^ /ef/™) m 



0.539 

0.72 

O' 

1.225 

0.80 

sF 

2.104 

0.89 

QlII 

2.820 

0.88 


creases with the incident irradiation energy (3.6%, 3.4%, 
3.1%, at 300, 600 and 900 eV, incident energies respec¬ 
tively). As extrapolated to larger energies, this fraction 
would decrease to 2.2% at 10 keV and 1.5% at 100 keV. 
The knowledge of the fraction of each type of defect with 
respect to the incident energy, as presented here, can now 
be used to estimate the total number of defects formed, 
based on the measurement on the concentration of oxy¬ 
gen vacancies only. 

It must be noted that the point defects that are formed 
can cause substantial modifications to the Si-0 network 
in quartz. For example, tri-cluster oxygen species are 
usually associated with edge-sharing Si tetrahedra, and 
three-member Si-0 rings are created by the relaxation 
caused by vacancies. These network defects constitute 
a major portion of radiation damage in quartz. These 
aspects, together with more complex structural damage 
in radiated quartz, would warrant further studies. 

B. Threshold displacement energy 

The threshold displacement energy is the energy 
needed to permanently displace an atom from its orig¬ 
inal lattice site. It is a critical input for the Monte Carlo 
method, which can be used to predict high-energy radia¬ 
tion damage for system sizes which are far too large to be 
handled by MD. This quantity is closely related to the 
atom types, the nature of bond and the crystal struc¬ 
ture, and is difficult to determine experimentally. We 
utilize the MD technique to evaluate threshold displace¬ 
ment energies for both silicon and oxygen atoms. They 
are obtained by accelerating a given atom in a quartz 
supercell with increasing energy until a permanent dis¬ 
placement is observed. The simulation is performed along 
different directions, at OK, and the results are presented 
in Fig. and Fig. It should be noted however, that 
the temperature should also have an effect on the thresh¬ 
old displacement energy. In particular, we keep in mind 
that radiation exposure will induce local heating, which 
could affect the stability of the atoms and, thereby, their 
threshold displacement energy. Such effect is beyond the 
scope of the present study. 

By plotting the threshold displacement energies for 
representative crystal planes, we observe a strong ori¬ 
entation dependence. In addition, the effect of the crys¬ 
tal symmetry is evident. Indeed, the minimum displace¬ 
ment energy is found to be 20 eV for oxygen along the 
<001> direction, and 30 eV for silicon along the <-110> 


direction among typical low-index orientations. Higher 
temperatures are likely to reduce the orientation depen¬ 
dence, resulting in a less pronounced orientation depen¬ 
dency. To predict the radiation damage in real applica¬ 
tions, where the material is usually randomly orientated 
with respect to the incident radiation, threshold energies 
are averaged over 20 uniformly spaced spherical angles. 
We find the average displacement energies of 28.9 eV for 
oxygen and 70.5 eV for silicon. The threshold displace¬ 
ment energy for oxygen has a similar value to that found 
in fused silica by previous simulation studies using a dif¬ 
ferent set of potentials p5l 15^ . This is expected as the 
local atomic structure is similar in both materials. Nev¬ 
ertheless, the fact that results obtained from different po¬ 
tential sets are consistent with each other demonstrates 
the robustness of the present method and suggests that 
it can generally be applied to study radiation-induced 
damage in silicate minerals or glasses. The value, how¬ 
ever, is considerably lower than the experimental value 
determined from optical measurements im, which may 
be related to the complex oxygen vacancy formation and 
annihilation process in irradiated quartz, as discussed in 
the previous section. We also note that the threshold en¬ 
ergy for silicon determined from our simulations is around 
twice as large as that of oxygen, agreeing with the phys¬ 
ical picture that twice as many bonds need to be broken 
to displace a silicon atom in comparison with an oxy¬ 
gen atom. Taken together, the validity of the present 
results, and their general nature suggests that MD sim¬ 
ulation can indeed provide a realistic atomic picture of 
radiation-induced damage in quartz. 


IV. CONCLUSION 

Radiation damage in quartz is studied by detailed MD 
simulations. By analyzing hundreds of cascade trajec¬ 
tories, we find that the ballistic events cause a range 
of structural damage. At the end of the cascade, over¬ 
coordinated Si and O species are found to be the main 
types of point defects, and their abundance strongly de¬ 
pends on the incident energy. Some under-coordinated 
defects also exist around the knock-on atoms, but their 
quantities remain limited, due to recombination and re¬ 
placement collisions. Oxygen vacancy defects are found 
to exist in a relatively low concentration, as compared 
with other defect types. In addition, this concentration 
decreases with the incident radiation energy. This sug¬ 
gests that the measurements of damage level based on 
E' center-related optical properties alone may not be a 
rigorous indicator of radiation damage. Structural dam¬ 
age also involves, in a considerable way, some modifica¬ 
tions of the Si-0 network, thereby creating small Si-0 
rings and edge-sharing Si tetrahedra. Such modifications 
could be the cause of radiation-induced amorphization, 
and should be investigated in further detail. We quantify 
the threshold displacement energy for both oxygen and 
silicon atoms from MD simulations. As expected, they 
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FIG. 6. Threshold displacement energy of oxygen along different directions in quartz crystal. 
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FIG. 7. Threshold displacement energy of silicon along different directions in quartz crystal. 


are found to be strongly dependent on the crystal ori¬ 
entation. The average threshold displacement energy is 
determined to be 28.9 and 70.5 eV for oxygen and silicon, 
respectively. These quantitative results enable improved 
assessments of radiation damage in quartz. 
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